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[l] The temporal stationarity of the teleconnection between 
the tropical Pacific Ocean and North America (NA) is 
analyzed in atmosphere-only, and coupled last-millennium, 
historical, and control runs from the Coupled Model 
Intercomparison Project Phase 5 data archive. The 
teleconnection, defined as the correlation between December- 
January-February (DJF) tropical Pacific sea surface 
temperatures (SSTs) and DJF 200 mb geopotential height, 
is found to be nonstationary on multidecadal timescales. 
There are significant changes in the spatial features of 
the teleconnection over NA in continuous 56-year 
segments of the last millennium and control simulations. 
Analysis of atmosphere-only simulations forced with 
observed SSTs indicates that atmospheric noise cannot 
account for the temporal variability of the teleconnection, 
which instead is likely explained by the strength of, and 
multidecadal changes in, tropical Pacific Ocean variability. 
These results have implications for teleconnection-based 
analyses of model fidelity in simulating precipitation, as 
well as any reconstruction and forecasting efforts that 
assume stationarity of the observed teleconnection. 
Citation: Coats, S., .1. E. Smerdon, B. I. Cook, and R. Seager 
(2013), Stationarity of the tropical pacific teleconnection to 
North America in CMIP5/PMIP3 model simulations, Geophys. 
Res. Lett., 40, 4927-4932, doi:10.1002/grl.50938. 

1. Introduction 

[ 2 ] The tropical Pacific Ocean impacts regional hydroclimate 
variability in the extratropics by means of wave propagation 
from areas of persistent precipitation and divergence anom- 
alies that are in turn forced by sea surface temperature (SST) 
variations [e.g., Sarachik and Cane, 2010], The preferred 
circulation responses to tropical Pacific SST forcing are called 
atmospheric teleconnections and depend on a Rossby wave 
response (see Trenberth et al. [1998], for a review) and 
subsequent interaction between the mean flow anomaly and 
transient eddies [Hoerling and Ting, 1994; Seager et al, 
2003, 2010; Harnik et al., 2010]. A quintessential feature of 
the El Nino-Southern Oscillation (ENSO) teleconnection is 
the Pacific jet stream shift over the western coast of North 
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America (NA), which a wealth of research has implicated as 
an important forcing of hydroclimate variability in NA on 
interannual to decadal timescales [e.g., Dettinger et al., 1998; 
Schubert et al., 2004a, 2004b; Seager et al., 2005a, 2005b, 
2008; Herweijer et al., 2006]. 

[ 3 ] The basis of our understanding of atmospheric 
teleconnections, between ENSO and NA or otherwise, is 
derived from the temporally limited observational record 
and associated reanalysis products [ Trenberth et al., 1998]. 
The stationarity of teleconnections is important for under- 
standing the impact of internal variability and forced 
responses within the climate system, an understanding that 
is fundamental to robust projections of regional climate 
change in a warming world. Within the approximately 100-year 
observational record, there is evidence that teleconnection 
patterns can vary considerably in space and time [e.g., Hu 
and Feng, 2001; Rajagopalan et al., 2000; Cole and 
Cook, 1998; Gershunov and Barnett, 1998]. A characteriza- 
tion of the stationarity of teleconnection patterns on 
decadal-to-centennial timescales, however, is not possible 
using only reanalysis and observational data because of 
their limited temporal extent. These limitations necessitate 
the use of alternative approaches to further characterize 
and understand teleconnection stationarity spanning decades 
and centuries. Toward such ends, multiple model products 
from the recent and ongoing release of the fifth Coupled 
Model Intercomparison Project (CM1P5) [Taylor et al., 2012] 
and third Paleoclimate Model Intercomparison Project 
(PMIP3) provide a state-of-the-art ensemble of model simu- 
lations for evaluating and testing teleconnection stationarity 
over multiple timescales. 

[ 4 ] Here we analyze three ensembles of fully-coupled 
Atmosphere-ocean General Circulation Models from the 
CMIP5/PMIP3 archive and compare them to the National 
Centers for Environmental Prediction/National Center for 
Atmospheric Research (NCEP/NCAR) reanalysis project 1 
[Kalnay et al., 1996] to assess the stationarity of the 
modeled ENSO-NA teleconnection. Our analyses suggest 
considerable variability in the simulated character of the 
teleconnection across the model ensemble on multidecadal 
timescales. These findings are crucial for assessments of 
model fidelity that often use regional teleconnection strength 
and spatial features as an evaluation metric, particularly 
with regard to simulated precipitation variability [e.g., Jospeh 
and Nigam, 2006; Smith and Chandler, 2010; Weare, 
2012; Langenbrunner and Neelin, 2013; Polade et al., 
2013], Furthermore, if realistic, these findings have impor- 
tant implications for forecasting [e.g., Barnston et al., 
2010] and reconstruction [e.g., MacDonald and Case, 
2005; D ’Arrigo et al., 2005] efforts that rely on the station- 
arity of teleconnection characteristics as defined by the 
observational interval. 
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Figure 1 . Teleconnection stationarity, as measured by the 
CPCS over NA, using teleconnection patterns estimated from 
the NCEP-NCAR reanalysis (see Figure 2) and the AMIP 
simulations. The first panel on the left is the CPCS value 
for each AMIP simulation. The second panel from the left 
plots the CPCS range for all 32 AMIP simulations. The third 
panel is the CPCS range for the 1979-2005 C.E. period in the 
1 6 coupled historical runs. The last panel on the right is the 
CPCS range for the reanalysis bootstrap experiment. Box 
plots represent the 25th and 75th percentiles of the data with 
the mean as the central line and the whiskers showing the full 
data range excluding outliers; outliers are marked with a 
cross (larger than 75th +1.5(75th — 25th) or smaller than 
25th — 1.5(75th —25th)). The number of ensemble members 
from each AMIP model are the values inset in the main panel. 

2. Data 

[ 5 ] All model output is from the CMIP5/PMIP3 archive. 
The three fully-coupled model ensembles are ( 1 ) 16 historical 
runs (1850-2005 Common Era (C.E.)) that have been forced 
with a common set of time- varying estimates of incident solar 
radiation, volcanic activity, and atmospheric constituents; 
(2) sixteen 500-year control simulations with constant 
preindustrial (pre-1850 C.E.) boundary conditions using the 
same models as in the collection of historical runs [Taylor 
et ah, 2012]; and (3) seven last-millennium (LM) simulations 
(850-1850 C.E.) forced with reconstructed time- varying 
exogenous forcing conditions [Schmidt et ah, 2011]. As a 
compliment to the fully-coupled runs, we use an uncoupled 
atmosphere-only ensemble comprising 10 sets of atmosphere 
model simulations with ensemble sizes between 1 and 10 from 
the Atmospheric Model Intercomparison Project (AMIP — a 
total of 32 simulations). Each ensemble member consists of 
an atmospheric simulation forced with the same observed 
SST boundary conditions for the period 1979-2008 C.E. All 
output from the above collection of simulations has been 
regridded to a common 2.5° x 2.5° latitude-longitude grid to 
match the reanalysis data set (for 66% of the simulations this 
represents a spatial degradation of the data). 

[6] The reanalysis data are taken from the National Oceanic 
and Atmospheric Administration NCEP/NCAR Climate Data 
Assimilations System 1 project [Kalnav et al., 1996], The 
monthly-resolved data span the period 1949-2012 C.E. on 
an equal area 2.5° x 2.5° latitude- longitude grid. The 
reanalysis product uses a frozen analysis/forecast system with 
data assimilation of past observations and is an established 
data set for assessment of interannual upper-air variability. 


3. Twentieth Century Teleconnection Variability 

[ 7 ] The analyzed model fields are surface temperature 
and 200 mb geopotential height. The latter was chosen over 
precipitation because it comprises a more spatially and tempo- 
rally homogeneous representation of the NA teleconnection, 
and because it is the ultimate driver of the precipitation vari- 
ability. For both fields, the climatology (calculated over the 
full simulations and the full length of the reanalysis) has been 
removed and the December-January-February (DJF) anoma- 
lies averaged. Winter averages were chosen because winter 
is the dominant period of precipitation forcing by ENSO in 
Western NA [Trenberth et al., 1998], 

[s] The teleconnection is defined by the grid point correlation 
between the DJF average 200 mb geopotential height fields and 
the DJF average Nino 3.4 (170°W-120°W, 5°N-5°S) SST index. 
The resulting correlation field indicates both the teleconnection 
strength and its spatial features. As an estimate of the observed 
teleconnection pattern, the correlation field was calculated for 
two periods in the NCEP-NCAR reanalysis: 1949-2005 C.E. 
(the overlapping period with the historical runs) and 
1979-2005 C.E. (the overlapping period with the AMIP runs). 
To assess the modeled teleconnection over NA, the centered 
pattern correlation statistic (hereinafter CPCS — following 
Santer et al. [1995]) was calculated between these reanalysis 
teleconnection patterns and the teleconnection patterns from 
the models over the NA region (160°W-50°W, 70°N-20°N). 
The CPCS is equivalent to a Pearson’s product-moment linear 
correlation between spatial values [Santer et al., 1995]. 

[ 9 ] Figure 1 characterizes the teleconnection variability 
across the collection of AMIP models and within the individual 
model ensembles. For each AMIP ensemble member, a 
teleconnection pattern was calculated and compared to the 
27-year period from the reanalysis using the CPCS. Because 
each simulation is forced with the same observed SSTs, the 
range in CPCS within each AMIP model ensemble represents 
an estimate of the impact of internal atmospheric variability 
on the circulation over NA. The calculated range was subse- 
quently compared to the pattern correlation range from the 
same period in the historical runs, which have coupled SSTs, 
and thus different SST variability and patterns. When all of 
the AMIP simulations are considered (32 in total), the 
CPCS range (excluding outliers) is about 0.35 (Figure 1). 
This is less than half the range over the same period 
(1979-2005 C.E.) in the coupled historical runs with 
dynamically evolving SSTs (Figure 1). Atmospheric noise, 
therefore, cannot be considered the dominant driver of differ- 
ences between modeled teleconnections overNA in the collec- 
tion of model simulations, and the relative differences in the 
range of the CPCS statistic are interpreted as arising from 
the SST-induced variations on the teleconnection over NA in 
the coupled simulations. For a further discussion of this analy- 
sis and interpretation, see the supporting information. 

[ 10 ] To evaluate the temporal stationarity of the observed 
teleconnection in the reanalysis data over a 27-year window, 
a bootstrap resampling of CPCS values in continuous 27-year 
segments was computed over the 56-year reanalysis record. 
The range in CPCS for 1000 of these resampled segments 
against the 1979-2005 C. E. reanalysis target is also plotted 
in Figure 1. The narrow range indicates that the character of 
the reanalysis teleconnection in the 27-year overlapping period 
between the AMIP and historical simulations is well repre- 
sented relative to the full 56-year period of the reanalysis. 
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Figure 2. Teleconnection stationarity, as measured by the CPCS over NA, using teleconnection patterns estimated from the 
NCEP-NCAR reanalysis and nonoverlapping 56-year segments from the 500-year control (solid boxes) and 1000-year forced 
LM runs (dashed boxes). Box plots indicate the 75th and 25th percentile of the CPCS statistic across the segments in the 
respective coupled model runs with the mean as the central line and the whiskers showing the full data range excluding outliers 
(larger than 75th +1.5(75th —25th) or smaller than 25th — 1.5(75th —25th)). The right four panels are the teleconnection 
pattern over 160°W-50°W, 70°N-20°N for the (first panel) reanalysis, (second and fourth panels) the most and least 
realistic segments, respectively, and the (third panel) model ensemble average. The colorbar range is -1-blue to l-red. LM 
ranges have been included as per their availability in the PMIP3 archive. Only models with 500-year control simulations were 
included, as a consequence, the GISS-E2-R and F-Goals models only have a LM simulation. 


4. Multidecadal Teleconnection Variability 

[l l] To investigate the nature of teleconnection stationarity 
over longer timescales, the control and LM runs were divided 
into 56-year segments to match the length of the NCEP/NCAR 
reanalysis record, generating a time-slice ensemble of 17 
members for each LM run and eight members for each control 
run. For each segment, the correlation between the Nino 3.4 
index and the 200 mb geopotential height field was calculated 
and compared to the 56-year (1949-2005 C.E.) pattern from 
the reanalysis, again using the CPCS. The range in the CPCS 
is thus interpreted as a measure of the temporal stationarity of 
the teleconnection within a given model. 

[ 12 ] The range of CPCS for the 16 control and seven LM 
runs is shown in Figure 2 and indicates a wide range in the 
teleconnection character within and between models. The 
models that have a stationary teleconnection (small range in 
CPCS values) in the control runs possess a similarly stationary 
teleconnection in the corresponding LM runs (e.g., Community 
Climate System Model version 4 (CCSM4)). The converse 
is also true (e.g., Model for Interdisciplinary Research on 
Climate Earth System Model (MIROC-ESM) and Max Planck 
Institute for Meteorology Earth System Model (MPI-ESM-P)). 
Transient forcing characteristics therefore do not appear to 
significantly impact the simulated teleconnection stationarity. 

[ 13 ] There is considerable spread in the CPCS range between 
models. For instance, the CCSM4 model simulates a stationary 
teleconnection that is consistently comparable to the reanalysis 
data. Canadian Earth System Model (CanESM) is likewise sta- 
tionary but less consistent with the reanalysis pattern, while the 
teleconnection simulated by the Commonwealth Scientific and 
Industrial Research Organization (CSIRO)-Mk3-6-0 model is 
neither stationary nor consistent with the reanalysis pattern. 
These observations are not explained by model resolution: 
while the CCSM4 model has the highest resolution of the 


simulations (0.9° x 1 .3°), CanESM is relatively low resolution 
(2.8° x 2.8°) and CSlRO-Mk3-6-0 is in between (1.9° x 1.9°). 
We therefore investigate below the potential dynamical links 
between variability in tropical Pacific SSTs and changing 
teleconnection characteristics over NA. 

5. Potential Dynamical Influences on Temporal 
Teleconnection Variability 

[ 14 ] What might cause the teleconnection over NA to be 
nonstationary? The dynamics of teleconnection variability 
within the models are not fully explored herein, but our anal- 
ysis suggests that a significant role is played by both the 
strength and spatial features of the SST anomalies in the 
tropical Pacific. 

[ 15 ] We use the Center of Pleat Index (CHI) \Giese and 
Ray, 2011] as a measure of the strength and location of 
SST anomalies in the tropical Pacific. The CHI statistic is 
analogous to the first moment of the SST anomaly field and 
provides both an amplitude and mean longitude for each 
ENSO event. Figures 3a and 3b plot average CHI amplitudes 
for each 56-year segment from the 16 control runs against the 
CPCS overNA during the same segment to assess the impact 
of the magnitude of ENSO events on the teleconnection. 
Segments with larger El Nino and La Nina events have 
higher values of the CPCS, although the connection is weak 
(Figures 3a and 3b). The same analysis was completed for 
CHI longitude and a much weaker relationship between 
eastward CHI longitude anomalies and high values of 
CPCS over NA was found (not shown). 

[ 16 ] Figure 3c investigates the impact of a model’s tropical 
Pacific SST variability on teleconnection stationarity by 
plotting the average variance of the Nino3.4 index for each 
time-slice ensemble against the range in the CPCS. There is 
a significant positive relationship between the ENSO amplitude 
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Figure 3. Average CHI amplitude for (a) La Ninas and (b) El Ninos in each 56-year segment plotted against the CPCS for that 
segment in the model control runs. The circle radius is proportional to the number of ENSO events in each segment, (c) The var- 
iance of the Nino3.4 index for the full control and LM simulations against the range in CPCS for that model simulation, (d) The 
range in the CPCS of the first EOF of tropical Pacific SSTs for the time-slice ensemble against the teleconnection CPCS range, (e) 
The CPCS of the first EOF of tropical Pacific SSTs from the top eight segments of each control run (see methods for ranking 
criteria) with that from the NCEP-NCAR reanalysis plotted against the CPCS of the NA teleconnection for that segment. The lines 
on each panel are the linear regression lines calculated using ordinary least squares. (See Table 1 for regression statistics.) 


and teleconnection stationarity, suggesting that large tropical 
Pacific SST variability allows the associated teleconnection to 
emerge above noise from secondary patterns of variability. 
For a discussion of this result in relation to the role of internal 
atmospheric variability in teleconnection stationarity, see 
the supporting information. 

[ 17 ] To detennine the impact of the stationarity of the 
modeled ENSO spatial patterns on teleconnection station- 
arity, we compute the first empirical orthogonal function 
(EOF) of the DJF tropical Pacific SST anomaly (defined as 
120°E-60°W, 30°N-30°S) for each segment of the control 
and LM simulations and compare them, using the CPCS 
statistic, to the same field from the reanalysis. This EOF 
was determined to be representative of ENSO by a strong 
correlation (average Pearson’s correlation coefficient of 
0.96 for all control and LM segments) between its principal 
component time series and the Nino3.4 index. In Figure 3d, 
the range in the CPCS between the model and reanalysis 
EOFs is compared to teleconnection stationarity to determine 


if models with large multidecadal variability in the spatial 
expression of ENSO have nonstationary teleconnections. 
There is a weak positive relationship, indicating that tempo- 
ral changes in the spatial character of ENSO exert a limited 
control on the stationarity of the NA teleconnection. 

[is] To further characterize the spatial nature of modeled 
ENSO, the first EOF of the DJF tropical Pacific SST anomaly 
was calculated for each control run using a 56-year sliding 


Table 1 . Pearson’s Coefficient of Determination (J? 2 ), P Values (Pval), 
and Mean Squared Error (MSE) for the Linear Least Squares Regression 
Derived From the Variable Associations Shown in Figures 3a-3e 


Figure 3 

R 2 

P Val 

MSE 

(a) 

0.068 

0.003 

0.022 

(b) 

0.063 

0.004 

0.040 

(c) 

0.336 

0.003 

0.175 

(d) 

0.096 

0.157 

0.000 

(e) 

0.399 

0.000 

0.000 
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window (incremented by 1 year, thus including overlap) and 
compared, using the CPCS statistic, to the same field from the 
reanalysis. The top eight segments for each simulation, as de- 
termined by the largest variance explained by the first mode of 
tropical ocean variability, were subsequently chosen for anal- 
ysis. This method of segment identification was done in order 
to maximize the EN SO-forced extratropical atmospheric sig- 
nal, using the rationale that segments with a dominant first 
mode of tropical Pacific Ocean variability will contain less 
“noise” from secondary patterns in the tropical Pacific 
Ocean. Additionally, there is a 0.86 correlation between 
the variance of the Nino3.4 index and the percent variance 
explained by the first EOF across all models, thus the iden- 
tified segments will also have the largest ENSO amplitude. 
The results are plotted in Figure 3e and indicate that the 
modeled and observed teleconnection will compare well 
when both strong ENSO variability exists and when the 
model’s ENSO spatial features match those of the reanalysis. 

[ 19 ] The statistical significance of the relationships in 
Figure 3 is investigated using Pearson’s correlation coeffi- 
cients. Despite our use of this metric for a quantitative descrip- 
tion of the variable associations, there is no reason to assume 
that the functional form is linear or that the relationship can 
be fully captured given the small sample size. While there is 
clearly a relationship between temporal and spatial changes 
in ENSO and the NA teleconnection, the correlations are 
weak (Table 1). This is not surprising given other sources of 
atmosphere-ocean variability. Nevertheless, the correlation in 
Figure 3e, our best attempt to isolate the signal of the ENSO 
spatial pattern and amplitude, indicates a moderate and statis- 
tically significant relationship. 

6. Discussion and Conclusion 

[20] The similarity of the ENSO-NA teleconnection 
between coupled model simulations and the observational 
record over the historical interval is an often used metric of 
model fidelity in reproducing coupled ocean-atmosphere 
dynamics [e.g., Joseph and Nigam, 2006; Smith and 
Chandler, 2010; Polade etal., 2013]. We nevertheless have 
shown that many models exhibit considerable variability in 
their teleconnection strength and character during different 
56-year windows of continuous simulations. These results 
suggest that analyses of teleconnection fidelity should be 
limited to atmosphere-only simulations forced with observed 
SSTs [e.g., Langenbrunner and Neelin, 2013]. Furthermore, 
considerable variability in modeled ENSO on multidecadal- 
to-centennial timescales has been reported in the literature 
[e.g., Wittenberg, 2009; Karnauskas et al., 2012; Ogata 
et al, 2013]. If such modeled behavior is dynamically 
realistic and drives changes in the atmospheric teleconnections, 
the 56-year reanalysis period may not be sufficient for 
adequately representing the character of the teleconnection over 
long timescales. In particular, the impact of nonstationary 
teleconnections on regional precipitation variability must be 
better understood, particularly given the range of reconstruction 
[e.g., MacDonald and Case 2005 and D ’Arrigo et al., 2005] 
and forecasting efforts (see Barnston et al. [2010] for a review 
of seasonal forecasting efforts) that rely on the stationarity of 
the observed teleconnection. 

[21] Our conclusions also have important implications for 
efforts to explain the hydroclimate history of western North 
America and its links to tropical Pacific climate variations. In 


particular, further characterizations of inter-model differences 
and outliers, as well as a more robust analysis of the coupled 
atmosphere ocean dynamics that drive temporal variations in 
simulated teleconnections, will provide improved understand- 
ing of teleconnection behavior and model performance. 
Furthermore, attempts must be made to relate teleconnection 
fidelity and stationarity to hydroclimate impacts on decadal 
to multidecadal timescales [e.g., Coats et al., 2013; also in 
preparation, 2013], Such efforts will improve characterizations 
of uncertainties in future projections of regional hydroclimate 
tied to large-scale teleconnections. 
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